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The ability of semi-classical transport models to correctly simulate Pauli 
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t^- , blocking and Bose enhancement is discussed. In the context of simple quan- 
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adequate to describe systems far from equilibrium. The relaxation of (1 + /) 
descriptions toward equilibrium is studied both in the context of simple mod- 
els, including a hot pion gas confined to a box. We give simple estimates for 
the characteristic number of collisions for a system to relax toward equilib 
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turn mechanical systems, it is shown that using (1 ± /) enhancements is in- 



rmm. 



Semi-classical simulations form the backbone of the phenomenology for high-energy and 
intermediate-energy heavy-ion physics. The semi-classical nature of such descriptions is 
justified if Compton wave lengths are smaller than the size of the reaction volume, and if 
phase space filling factors, /, are much less than unity. The first condition is generally 
fulfilled whenever incident beam energies exceed ~ 40 • A MeV. The second condition is 
violated for lower energy collisions by nucleons, as excitation energies are less than the 
Fermi energy until beam energies exceed ~ 100 -A MeV. For relativistic heavy-ion collisions, 
E/A > 10 GeV, phase-space-filling factors for pions also approach unity. For these two 
regimes, quantum statistics must be incorporated into simulations. 
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Nordheim [|IJ first pointed out that modifying the Boltzmann transport equation by 
incorporating phase-space occupation factors (1 ± /) into the collision integral, showing the 
phase space distribution functions have the Bose-Einstein or Fermi-Dirac distributions as 
their respective equilibrium solutions. The usual way to accomplish this is by modifying the 
cross sections from their in- vacuum values (y v ^L, cd to the in- medium value |@-§J: 

a™£ cd = <>l£> cd .{l±f c ){l±U). (1) 

Given a system with measured occupations, a scattering will choose its final state ac- 
cording to the weights described above. This problem was originally cast in terms of the 
Einstein A nm and B nm coefficients ||, which is the correct treatment given measurements 
are made after each emission. However, when measurements are not made between indi- 
vidual scatterings; correct treatment of Bose and Fermi statistics requires the incorporation 
of many-body effects. Our principal objective here is to compare models where full n-body 
symmetrization has been accounted for with the (1 ± /) method. 

We start by examining the Einstein problem: If in free space the rate at which quanta 
are scattered into a cell k is a k and the rate at which quanta are scattered out of the cell 
is Pkfk, the rate equations for the phase-space occupation f k in and out of the medium are 
given by: 



JMvac 

\ = a* - /r%- (2) 

J-Cmed 

h -a k (l±fr d )-fr d Pk. (3) 



dt 

This results in equilibrium values of 

nr = f (4) 

rmed Pk /r\ 

Jk - 1 3- Ofc ■ K°) 

If the free space values of a k and f3 k result in a Boltzmann factor, a^/Pk — exp(— {E^ — 

n)/T), the in- medium result is a Bose-Einstein or Fermi-Dirac form. Thus one can be 
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confident that a simulation that incorporates (1 ± /) corrections will ultimately obtain 
correct equilibrium behavior given sufficient time such that gain and loss terms balance. 

In a relativistic heavy-ion collision (E/A ~ 200 GeV) the reaction can be viewed as 
having two stages. During the creation stage pion production takes place, while during 
the reinteraction stage expansion and cooling occur, during which the pion number stays 
roughly constant. We divide our model into two stages correspondingly. For the first stage 
we consider only production. Given the large energies of the colliding nucleons, reformation 
of the energetic initial states is negligible and one can consider only formation processes. In 
the absence of Bose effects, we assume the probability of having created n pions is Poissonian: 

P fc V,0=exp(- % (i))^. (6) 

71/. 

Fully including symmetrization, the probability of creating n pions is enhanced by a factor 
of n\ assuming the time t is after emission has stopped and that no measurements were made 
at intermediate times. This enhancement can be considered as either the contribution of 
interference terms or as the amplitude of a matrix element with n identical Boson creation 
operators, (0\a n aJ n \0). 

Pt ull (n,t) = ±P° k (n,t)n\ (7) 
= (1 - Vk(t))Vk(t) n 

where Z k is a normalization. 

If the state is measured at each intermediate time t during the emission, (1 ±/) methods 
are justified. The probability of seeing n particles at time t, Pj, 1+ ^ (n,t) is then determined 
by the equation: 

j t Pr\n, t) = *tj$. {„ ■ Pi' + '\n,t) - (1 + n) ■ Pt f, (n + 1, *)} (8) 
with the solution 

Pi 1+f \n,t) = (l- lk (t)) lk (tr (9) 
7 fc (*) = 1 - exp(-77 fc (t)) 

3 



Simulating systems populated according to (1 + /) weights results in dendritic growth 
of phase-space occupation. The first few particles can bias the behavior of all the following 
particles, leading to large fluctuations in the populations. Whereas, when no measurement 
is made until the end of the emission, all of the emissions are affected equally by the sym- 
metrization. The average occupation is lower for the (1 + /) case, Eq. (^j), than for the case 
where full n-body symmetrization weighting is applied at the end of the emission as in Eq. 
(0), but certainly higher than for the case where symmetrization is neglected, 



(n°k) = Vk{t=oo) 
lk(t = cx 

,fuih Vk(t=oo) 
\ n k I 



= S (10) 



1 - Vk(t = oo) 

If the function rjk has a thermal form, % = exp(— (E^ — u)/T), the full n-body result 
corresponds to the thermal Bose-Einstein form. The additional created particles account for 
the difference of the Bose-Einstein form as compared to the Boltzmann form for the spectra. 
For the (1 + /) case, the population of the different levels of energy Ek can not be described 
by one temperature and chemical potential. 

Fig. [TJ(a) illustrates the three solutions for the mean occupancy. There is more than just 
a quantitative difference. The (1 + /) result never diverges as compared to the full n-body 
case which results in an infinite occupancy when rjk exceeds unity. The (1 + /) result even 
differs to the second power of 

One can solve the same models for fermions. In this case (1 + /) terms are replaced by 
(1 — /) terms in Eq. (|8|) and the n\ in Eq. ([?]) is replace by zero for n > 1. The results look 
very similar for the mean occupations. 

<4 W) > = 7k(t = oo) (11) 

/full, _ Vk(t = oo) 
\ n k I 



1 + Vk(t = oo) 

The difference is illustrated in Fig. [TJ(b) which shows that (1 — /) models can modestly 
underestimate Pauli-blocking effects. 



Even if (1 + /) models grossly underestimate initial populations of low-energy states, 
rescattering effects can lead to equilibrium. If sufficient scattering occurs to produce equi- 
librium, the only fault of the (1 + /) models will lie in the prediction of the overall multi- 
plicity. For this reason we wish to investigate the amount of scattering required to produce 
equilibrium. 

We consider a closed system with jVj levels labeled by k where each level has a degeneracy 
d k and base weights r] k . For initial conditions we fill the levels with N particles according 
to the (1 + /) prescription, neglecting any loss terms. Collisions are modeled by randomly 
removing one of the N particles and placing it into a level according to the weight df.(l + 
n k )r] k . After sufficient collisions populations approach equilibrium. 

The equilibrium population for fixed N systems can be calculated by using diagram- 
matic methods ||. Although these method have only been used to study multi-particle 
symmetrization of outgoing states, they can also be used to calculate equilibrium popula- 
tions. Ref. || describes how to calculate the probability of emitting iV particles in states 
k\---k N given the source function S(p,x), where S(p,x) is the probability of emitting a 
particle from space-time point x into final state p neglecting symmetrization. Under the 
assumption that the momentum dependence in S is not correlated to the space-time depen- 
dence, S(p,x) = S(x)rj{p) the expression in Ref. |§ becomes: 



P(ki,---k N ) = — / dx 1 ---dx N rj kl ---r] kN S(x 1 )---S( y x N )\U(x 1 ,---x n ;k 1 ,---k N )\ 2 , (12) 



where U is the evolution matrix from the creation of pions at x to their measurement in 
states k. If we assume that the system is finite with discreet eigenstates, and if we assume 
that there is no x dependence aside from confining the region of emission, the expression 
above becomes: 



The squared quantity is the squared n-particle wave function with N\ terms. The same 
equation applies for an open (non-bounded) system if the emission occurs all at one time, 




P(k h ---k N ) = - J 



dxi ■ ■ -dx N i] kl --- T]k N \(xi, ■ ■ -x n \ki, ■ ■ ■ 




(13) 
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except that the labels refer to continuum states rather than discreet eigenstates. By making 
the substitution rji = exp(—Ei/T) one sees that Eq. ([13]) is equivalent to a thermal trace, 
hence it is the equilibrium answer. 

Methods for accounting for all N\ terms in the outgoing wave function were devised 
in Ref. ||. Translating the method of Ref. |J to the case of discreet states results in the 
prescription below for calculating (rik) given rjk and N. First define quantities G n (k), C n 
and uj{n). 

G n (k) = (il)"- 1 ^ (14) 
C n {k) = -Y,G n {k) 



/"til riii ri 

f(n) = 2^ [[ 



il— In 1 z n 

The sum over i\ ■ • ■ is constrained such that the overall order is N, i\ + 2%2 + 3«3 + ■ • ■ = 
N. An algorithm for calculating the sum in the expression of u(n) is also given in Ref. |§. 
The expression for the equilibrium occupations is: 

. G 1 (k)u(n-l) + G 2 (k)u(n-2) + ---G n ^(k)uj(l) + G n {k) 

\nk)equil - U(W) ' ' ' 

Now we return to the evolution of (n^) averaged over many samplings. Fig. Q(a) shows 
the occupation as a function of the number of collisions for the case where there are 8 
levels with filling weights rji = exp(—k ■ AE/T), with AE = T/8, and 40 particles. Given 
the large number of particles in a limited phase space one can easily see that the initial 
conditions underpopulate the ground state, but that a sufficient number of collisions restores 
equilibrium. The dashed lines in Fig. 0(a) represent the equilibrium occupancies calculated 
above. In order to make a detailed study of the approach to equilibrium Fig. 0(b) shows 
the discrepancies from equilibrium as a function of the number of collisions on a logarithmic 
scale. This allows the extraction of characteristic scales for the exponential approach to 
equilibrium. 

One striking aspect of Fig. |2](b) is that as one approaches equilibrium the lines become 
parallel. The inverse slope of this line corresponds to the characteristic number of collisions 



for approaching equilibrium r d . In order to understand how the parameters (level density, 
temperature, ground-state population) determine r^, we linearize the evolution equation and 
treat fk as a continuous variable. This allows us to quickly find r d for a variety of systems, 
which will lend insight into what characteristics determine r^. Ignoring the discreet nature 
of f k , 

dfk (1 + fk)Vk fk 



dr Z N' 

Z = Nj2(l + fk)Vk. 



(16) 



The normalization factor Z which is a function of the distribution f k , ensures that particle 
number is conserved. For small deviations from equilibrium, one can expand Z keeping the 
lowest order terms in Sf k = fk — fk to obtain: 

dSfk Sf k , f k « 5f 



+ 



Y^Yfr- (17) 



dr 1 + f? N j 1 + fj 
One can solve the N coupled differential equations with matrix methods. There are 
N — 1 solutions given the constraint J2k$fk = 0. Each solution has the form: 

M««p{-^£}. (18) 

Each eigenvector has a different value of r d . We find that the slowest decaying eigenvector 
is the one where the ground state is furthest from equilibrium. For thermal systems we also 
find one can estimate r d m T/ (E 2 — E{) where E x and E 2 refer to the energies of the ground 
state and the first excited state. This approximation works to better than 20% accuracy 
whenever the first state has an occupation much greater than unity. This is true even 
when the second state is degenerate. Thus in the case of a nearly degenerate ground state, 
equilibration is very slow. 

The simulation of the eight-level model described above can be easily modified to handle 
more sophisticated systems. To that end we consider the eigenstates of a cube of length L = 
10 fm on a side with a base filling factor, rji = exp(—Ei/T), with T = 175 MeV. Assuming 
the eigenstates are determined by a relativistic dispersion: 



El M = m 2 + ^(kl + kl + kl) , (19) 

we simulate the behavior of (fi) for 200 particles. Fig. ||] shows symmetrization enhance- 
ments to spectra. The symmetrization enhancement is defined by the average phase-space 
filling for states of that energy divided by exp(—Ei/T). The lowest curve demonstrates the 
symmetrization enhancement after states have been filled according to the (1 + /) prescrip- 
tion with no collisions. The figure shows that after five collisions per particle symmetrization 
enhancements for small energies are half way to equilibrium values. The estimate of men- 
tioned above for this system is six collisions per particle. 

In heavy-ion collisions the number of collisions pions feel is on the order of a half dozen 
for sulphur on lead collisions at 200 GeV incident laboratory energy. For the upcoming 
lead beam that number could roughly double. Thus we see that for larger systems, (1 + /) 
theory should be able to recover from its shortcomings and predict a spectra in line with what 
would be expected had all n-body effects been included. However, the multiplicity would 
still be underestimated. The amount by which the multiplicity might be underestimated 
would be similar to the amount of emission due to stimulated emission. Fig. [I] shows that 
the extra multiplicity could possibly be very large if values of rjk become larger than unity. 
Although a real system can not emit an infinite number of particles which is possible with the 
Poissonian model used above, the possibility of a significant fraction of pions being due to 
enhanced emission is possible. Phase-space occupations should be even more overpopulated 
for conditions created by the relativistic heavy ion collider (RHIC), scheduled to come on 
line in 1999. 

In proton-induced collisions, conclusions are simpler. For these systems reinteraction of 
hadrons is not so important. Here (1 + /) models are completely inappropriate and full 
n-particle symmetrization effects must be incorporated. This issue has been addressed in 
the context of centauro phenomena f7HlQfl. 

The possibility of creating super-radiant pion sources is one of the most exciting prospects 
of relativistic-heavy-ion physics, cosmic-ray physics and particle physics. Symmetrization 
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phenomena should manifest themselves in a variety of observables: multiplicity distributions, 
isospin distributions and spectra. Making confident predictions of collective symmetrization 
phenomena is inherently difficult due to the unstable nature of enhanced emission. We 
have shown here that (1 + /) based simulations can grossly underpredict such phenom- 
ena. Should super-radiant conditions be attained in relativistic heavy-ion collisions, one of 
the most daunting theoretical challenges will be finding a way to incorporate many-body 
symmetrization into descriptions where emission and collisions are treated realistically. 
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FIGURES 

FIG. 1. a: Average population of a Bosonic state as a function of the Poissonian parameter rj 
for the cases of neglected symmetrization (dotted line), (1 + /)-enhancement (dashed line), and 
full n-body enhancement (solid line), b: Same comparison for a Fermi gas; (1 — /) suppression 
(dashed line), full n-body suppression (solid line). 

FIG. 2. Eight level system with 40 particles and level spacing AE = T/8. 
a: Average populations of the eight levels as a function of the number of collisions per particle. 



Dashed lines: exact equilibrium values computed from Eq. 15. Solid lines: (l + /)-type simulations, 
averaged over 3 x 10 6 events. 

b: Deviation of the (1 + /)-type simulations from equilibrium populations. After many collisions 
the lines become parallel. Statistical noise causes the ragged behavior for small deviations. 

FIG. 3. Populations for states, divided by exp(— Ei/T), as a function of the energy for a system 
of 200 relativistic particles confined to a 10 fm cube at a temperature of 175 MeV, obtained from 
(1 + /)-type simulations and averaged over 10 5 events. Initial populations (circles), after one 
(crosses), two (diamonds), five (squares), and fifty (upright crosses) collisions per particle. After 
many collisions populations have approached equilibrium values (dashed line). 
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